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Behavior of two-time autocorrelation during the phase separation in solid binary 
mixtures are studied via numerical solutions of the Cahn-Hilliard equation as well as 
Monte Carlo simulations of the Ising model. Results are analyzed via state-of-the-art 
methods, including the finite-size scaling technique. Full forms of the autocorrelation 
in space dimensions 2 and 3 are obtained empirically. The long time behavior are 
found to be power-law type, with exponents unexpectedly higher than the ones for 
the ferromagnetic ordering. Both Chan-Hilliard and Ising models provide results 
consistent with each other. 
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I. INTRODUCTION 


Properties of a nonequilibrium system change with growing age 1]. Understanding of such 
aging phenomena is of fundamental importance in all branches of science and technology. There 


have been serious activities on this issue concerning 


especially in problems related to domain growth 


lvmg 


, 3] as well as nonliving matters, 


14] and glassy dynamics 


15 


19]. Among 


other quantities, aging phenomena is studied via the two-time autocorrelation function 


C(t,t w ) = {il>(r,t)il>(r,t w )) - ('ip(r,t))('ip{r,t w )). 


4] 


( 1 ) 


In Eq. (ITj). ib is a space (r) and time dependent order parameter, t w is the waiting time or age 


of the system and t (> t w ) is the observation time. 

[ 20 ], 


In phase ordering systems [20|, though time translation invariance is broken, C(t,t w ) is 
expected to exhibit scaling with respect to t/t w . Important examples are ordering of spins 
in a ferromagnet, kinetics of phase separation in a binary (A + B) mixture, etc., having been 
quenched to a temperature (T) below the critical value (T c ), from a homogeneous configuration. 
Though full forms are unknown even for very simple models, asymptotically C(t, t w ) is expected 


to obey power-law scaling behavior as 


4, 6] 


C(t , t w ) ~ x A ; x = Ijt 


( 2 ) 


In Eq. (|2j), £ and £ w are the average sizes of domains, formed by spins or particles of similar 
type, at times t and t w , respectively. Typically £ and t are related to each other via power-laws. 


For nonconserved order-parameter dynamics, e.g., ordering in a ferromagnet, such sea’ 


has been observed and the values of the exponent A have been accurately estimated 
different space dimensions d. There the exponents follow the bounds 

d x 
2 


6, [l4] in 


mg 


( 3 ) 
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predicted by Fisher and Huse (FH) 4]. In kinetics of phase separation in solid mixtures, 
for which the order parameter is a conserved quantity, the state of understanding is far from 
satisfactory, due to theoretical as well as computational difficulties. There exist reports [11] of 
violation of scaling with respect to t/t w . The latter observation, our results indicate, is due 
to the fact that scaling is achieved for t w ^ 1. Access to such long time is constrained by 


inadequate computational resources. This difficulty, in some studies 


13], might have led to 


the conclusion about incorrect values o 


In an important work, Young et al. 


A. 


put a more general lower bound on A, valid irrespective 


of the conservation of the total order parameter, as 


A > 


/3 + d 


( 4 ) 


where f3 is the exponent for small wave-vector power-law enhancement of equal time structure 
factor which, depending upon the dynamics, becomes important for t w 1, as stated below. 
In nonconserved dynamics, (3 = 0 and so the lower bound in Eq. ([3]) is recovered. For conserved 


order parameter dynamics, on the other hand, (3 = 4 in both d = 2 and 3 at late time, 
the upper bound in Eq. ([3]) is violated. Simulations of the Cahn-Hilliard (CH) equation 


hus, 


20 ] 


d*l>(r, t) 
dt 


= -V 2 


ip(r, t ) + V 2, 0(r, t ) — ijj 3 (r, t ) 


( 5 ) 


by Young et al. 


7], observed A > 3 in d — 2, consistent with Eq. ([4]). From these simulations, 
the authors, however, did not accurately quantify A; scaling of C(t,t w ) with respect to t/t w was 
not demonstrated; focus was rather on the sensitivity of the aging dynamics to the correlations 
in the initial configurations. Situation is far worse in d — 3, with respect to the CH equation 


as well as the Ising model 


20 ] 


H = ~JJ2 SiS i' Si = ±1; J> °- 


<ij> 


( 6 ) 
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In this paper, we study both CH equation and the Ising model, used for understanding 
diffusive phase separation as in solid mixtures, in d — 2 (on regular square lattice) and d = 3 (on 
simple cubic lattice), via extensive simulations, to quantify the decay of C(t,t w ). We observe 
scaling of C(t,t w ) with respect to x in which the power-law of Eq. fl2]) is realized for large x. 


Via computations of the instantaneous exponent 21 23] 


A,; = 


d\n[C(t,t w )] 
d In x 


(7) 


and application of the finite-size scaling technique 24 . 


25], we find that A ~ 3.6 in d — 2 


and ~ 7.5 in d — 3. Though these numbers respect the bounds in Eq. (J3J), the high value 
in d = 3 is surprising. But this comes from both CH and Ising models, from various reliable 
analyses. Furthermore, a general analytical form for the full scaling functions has been obtained 
empirically. 


II. METHODS 


One numerically solves the CH equations on a regular lattice, usually via Euler discretization 


method. With the Ising model, phase separation 
via the Kawasaki exchange Monte Carlo (MC) 


nineties in a solid binary mixture is studied 


25] simulations, to be referred to as KIM. An 


up spin (Si = +1), for this problem, may correspond to an A particle and a down spin (—1) to a 
B particle. In this MC scheme, one randomly chooses a pair of nearest neighbor spins and tries 
;heir position exchange. The moves are accepted according to standard Metropolis algorithm 


251] . Due to the coarse-grained nature of the CH equation, as opposed to the atomistic Ising 


model, one can explore large effective length in simulations. The order parameter in Eq. 


corresponds to a coarse-graining 


26] of the Ising spins, typically over the equilibrium correlation 


length £. Then, a positive value of ijj means an A-rich region and for a H-rich region, -0 will 
have a negative number. For the calculation of C(t,t w ), we have used binary numbers +1 and 
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-1, for both the models. 

The average domain length, £, in our simulations was measured from the first moment of 


domain size distribution, P(£ d ,t), as 


23] 


£ = / i d P(£d,t) M d , 


( 8 ) 


where £ d is the distance between two successive domain boundaries in any direction. Throughout 
the paper, all lengths are presented in units of the lattice constant a. In MC simulations, 
time is counted in units of Monte Carlo steps (MCS), each MCS consisting of L d trial moves, 
where L is the linear dimension of a periodic square or cubic system. In CH equation, t is 


expressed dimensionless units 


27]. All results are presented after averaging over at least 50 


initial configurations, for quenches from random initial configurations to T = 0.6 T c . 


III. RESULTS 

In Fig. Oja), we present the plots of C(t,t w ), vs x, for different values of t w , from the 
solutions of CH model in d — 2. As seen, one needs large enough value of t w to observe 


appropriate scaling behavior, compared to ordering in ferromagnets 14]. Between the two data 


sets with largest values of t w , the deviation from each other, for large x, is due to the finite-size 
effects. Similar plots for the d — 3 CH model are presented in Fig. [0(b). Here all t w s are large 
enough, providing good scaling. Again, deviation from the master curve, starting at different 
values of x for different t w , are primarily related to the finite-size effects. In both these figures, 
Q](a) and [D]b), the system sizes are kept fixed, only the values of t w are varied. A similar 
observation, with respect to the above mentioned deviation for different choices of t w , can be 
made, when, for same value of t w , data are presented for different system sizes. 

In the scaling parts, both in Fig. Oja) and Fig. [0(b), continuous bending is observable, in 
these log-log plots. Thus, power-laws, if exist, carry corrections. The solid lines in these figures 
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FIG. 1. (a) Autocorrelation function, C(t,t w ), from the d = 2 Cahn-Hilliard model, are plotted vs 
x (= £/£ w ), for different values of t w . The solid line there corresponds to a power-law decay with 
exponent 3. (b) Same as (a) but for the d = 3 CH model. The solid line there has a power-law decay 
exponent 3.5. The system sizes used are L = 256 (d = 2) and 200 (d = 3). 

are power-law decays with exponents 3 and 3.5, respectively, corresponding to the bounds in 
Eq. ((4j). For large x, simulation data in d — 2 appear reasonably consistent with the bound. 
The asymptotic exponent, in d — 3, on the other hand, appear much higher than 3.5. 

With the expectation that power laws indeed exist, in Fig. [2] we present plots of instanta¬ 


neous exponents 


, |l4j, 21 23] A i, for both d = 2 and 3, vs 1/x. In addition to providing A, 


from the extrapolations to x = oo, such exercise may be useful for obtaining crucial information 
on the full forms of C(t,t w ). For d = 2, the data are obtained for t w — 5 x 10 3 , and for d — 3, 
the data correspond to t w = 10 3 . In both the cases, the results appear reasonably linear p, [14]. 
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FIG. 2. Instantaneous exponents A* are plotted vs 1/x. Results are shown only from the solutions of 
the CH equations, in both d = 2 and 3. The solid lines are guides to the eye. The d = 2 data are for 
t w = 5 x 10 3 with L = 400. In d = 3 the numbers are 10 3 and 200. 


The solid lines there are extrapolations to x = oo, accepting the linear trends. These indicate 
A ~ 3.60 in d — 2 and ~ 7.80 in d — 3. Again, while the value in d — 2 is consistent and 


close to the bound of Yeung 


7] et ah, the observation of surprisingly high number in d — 3 


is certainly interesting. We intend to obtain more accurate values via appropriate finite-size 


scaling analyses 


24, 


25]. This is considering the fact that the choice of the regions in Fig. [2j 


for performing least-square fitting, is not unambiguous due to finite-size effects and strong sta¬ 
tistical fluctuations at large x. Also, for very small x (data excluded), there is rapid decay of 
C(t,t w ) related to the fast equilibration of domain magnetization. 


Since the corrections to the asymptotic decay laws are seen to be strong for finite x, a 
reasonable idea about the full forms of the decays is essential for accurate finite-size scaling 
analyses. Those, however, are nonexistent in the literature. Here we obtain the forms empiri¬ 
cally. Assuming power-law behavior of the data sets in Fig. [2l we write 


Ai = A - —, 

x Y 


( 9 ) 












1.1 



FIG. 3. Finite-size scaling plot of C(t,t w ) from d = 2 CH model. The scaling function Y is plotted 
vs y , using data from different system sizes. The optimum collapse of data, the presented one, was 
obtained for A = 3.47. The inset presents the same exercise for the CH model in d = 3. Here the 
value of A is 7.30. See text for values of t w . 

where A c and 7 are constants. Combining Eq. (0 with Eq. (|7|), we obtain 

C(t,t w ) = C 0 exp ( - ( 10 ) 

Co being a constant. For finite-size scaling analysis, one needs to introduce a scaling function 

Y(y) = C(t,t w )ex p (-^)x A ; y = L/t. (11) 

For appropriate choices of A c , 7 and A, one should obtain a master curve for Y, when data from 
different system sizes are used. The behavior of Y should be flat in the finite-size unaffected 
region and a deviation from it will mark the onset of finite-size effects. 

By examining the data in Fig. 0] (also see Fig. 01(b) for KIM), we fix 7 to 1. In the main 
frame of Fig. 01 we show a finite-size scaling plot for data from the d = 2 CH model. The 
presented results correspond to best collapse, obtained for A c = 2.25 and A = 3.47. The value 
of t w used for all the data sets is 10 4 . A similar exercise for the d = 3 CH data is presented 
in the inset of Fig. 0J Again the data collapse looks quite reasonable and was obtained for 
A c = 5.1 and A = 7.30. The value of t w , in this case, was set to 10 3 . The reason behind 
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choosing smaller value of t w in d — 3, than in d — 2, is computational difficulty. It is extremely 
difficult to accumulate data for further decades in time, starting from very high value of t w , 
particularly in d = 3. Nevertheless, this chosen value of t w falls within the scaling regime. Note 
that for similar temperatures, amplitude of growth is larger in d = 3 and the scaling of C(t, t w ) 
is related more closely to the value of l w . In this connection we mention that the longest run 
lengths (associated with largest systems) for the CH model are t — 2 x 10 6 and 2 x 10 5 in d — 2 
and 3, respectively; for the Ising model these numbers are 5 x 10' and 4 x 10 6 . 



FIG. 4. (a) Same as Fig. Uta) but for the d = 2 Ising model, (b) Same as Fig. [2] but for Ising model. 
The oscillatory behavior of large x data in d = 2 is statistical fluctuation. In d = 2 the results are 
from L = 512 and in d = 3, we have used L = 100. 


We now move to present results from KIM. In Fig. IU(a) we show the autocorrelations from 
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different values of t w in d = 2, for L — 512. Scaling is poor for t w below 10 4 MCS and so those 
results are excluded. Despite strong statistical fluctuations, it is recognizable that the decay 
of C(t,t w ) in the latter part is on the higher side of the bound in Eq. (J3J), represented by the 

solid line. 



FIG. 5. C(t, t w ) are plotted vs x , for the CH model in d = 2 and 3. Inset shows corresponding results 
for the KIM. The solid lines are fits to the form in Eq. m, with 7=1. The t w values are mentioned 
on the figure. We have discarded data suffering from finite-size effects. The system sizes are L = 400 
and 512 for CH and Ising models in d = 2, whereas, L = 200 and 100 in d = 3. 

In Fig. 0](b) we show the instantaneous exponents for the Ising model in d — 2 and 3, 
vs 1/x. In each dimension, we have included two values of t w from the scaling regime. While 
results for different t w s, in a particular dimension, are consistent with each other, f ini te-size 
effects appear earlier for larger value of t w , as expected. Thus, for extrapolations to x = oo, 
data sets with smaller t w are used. This exercise provides A ~ 3.60 and ~ 7.30 in d = 2 and 
3, respectively. These values are in agreement with the ones obtained for the CH model via 
various methods of analysis. 

Finally, in Fig. [5] we show the fits of the simulation data to the form in Eq. (fTUlh The main 
frame is for the CH model in d — 2 and 3, whereas the inset contains similar results from the 
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KIM. The fits look quite satisfactory. This exercise provides A = 3.55 and 3.76 in d = 2 for CH 
and Ising models, respectively. The numbers in d = 3 are 7.64 and 7.37. 


IV. CONCLUSIONS 


In conclusion, we have studied aging dynamics for the phase separation in solid binary 
mixtures via Cahn-Hilliard and Ising models. Results for the two-time autocorrelation, C(t, t w ), 
are presented from simulations in both d — 2 and 3. Decays of C(t,t w ) appear power law in 
large x limit. The exponents for these power laws were obtained via various different analyses, 
including finite-size scaling, which is new for this purpose. For the finite-size scaling analysis, full 
forms of the autocorrelations were essential which we obtained empirically. All these methods 
provide consistent values of the decay exponent A for different models. These are A ~ 3.6 in 
d — 2 and A ~ 7.5 in d = 3, within 5% error. 


Very high value of A in d — 3, far above the lower bound in Eq. flJJ), can be due to the fact 
that domains are more mobile in this dimension than in d — 2. From the numbers obtained 
in d — 2 and 3, it may be tempting to predict a dimensionality dependence as A = f(d — 1) 
with / ~ 3.75. However, we caution the reader not to jump into such conclusion. Even though 
the influence of the dimension d appear more important than in the bound of Eq. (HI), the 
contribution of (3 is significant, particularly at lower dimension. Results from other dimensions 
are necessary to make such a conclusion. In d = 1 one should exercise the caution that /3 (= 2) 
has a different value 28]. 
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